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Abstract — We report progress in the development of a 
model-based hybrid probabilistic approach to an on-board 
IVHM for solid rocket boosters (SRBs) that can 
accommodate the abrupt changes of the model parameters in 
various nonlinear dynamical off-nominal regimes. The work 
is related to the ORION mission program. Specifically, a 
case breach fault for SRBs is considered that takes into 
account burning a hole through the rocket case, as well as 
ablation of the nozzle throat under the action of hot gas 
flow. A high-fidelity model (HFM) of the fault is developed 
in FLUENT in cylindrical symmetry. The results of the 
FLUENT simulations are shown to be in good agreement 
with quasi -stationary approximation and analytical solution 
of a system of one-dimensional partial differential equations 
(PDEs) for the gas flow in the combustion chamber and in 
the hole through the rocket case. 

The low-dimensional performance model (LDPM) of the 
fault is derived by integrating a set of one-dimensional 
PDEs along the axis of the rocket. The LDPM is used to 
build a model-based fault diagnostic and prognostic (FD&P) 
algorithm for the case breach fault. In particular, two 
algorithms are introduced. The first algorithm is based on 
the self-consistent algorithm that solves the LDPM in a 
quasi-adiabatic approximation, when the pressure and 
density follow adiabatically dynamics of the propellant 
burning, melting and burning of the metal case, ablation and 
erosion of nozzle and insulator. The second algorithm is 
based on the dynamical inference method [l]-[3] of the 
system of stochastic differential equations of the LDPM. 

The parameters of the HFM model and of the LDPM are 
tuned to reproduce the results of recent experiments of the 
rocket firing with the case breach fault in the forward 
closure. The FD&P is then applied to illustrate real-time 
diagnostics of the model parameters and prognostics of the 
SRB internal ballistics. All the algorithms discussed in this 
paper were verified using experimental data as will be 
discussed elsewhere. 

The accuracy of the algorithm and the possibility of its 
application to FD&P for other SRB fault modes are 
discussed. 
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Nomenclature 



p 


= 


gas density 


p 


= 


gas pressure 


T 


= 


gas temperature 


u 




gas velocity 


v m 




velocity of metal melting front 


v„ 




velocity of nozzle ablation front 


c 




sound velocity 


M 




Mach number, M = u/c 


c v 




specific heat for constant volume 


c p 




specific heat for the constant pressure 


7 




ration of specific heats y= Cp/c v 


I 




perimeter of propellant cross-section 


e, 




total energy of the gas 


h, 




total enthalpy of the gas 


h 




perimeter of hole cross-section 


r„ 




radius of nozzle throat 


s, 




cross section of the nozzle throat 


r h 




minimum radius of leak hole 


Sh 




cross-section of hole throat 


'met 




radius of hole in metal case 


r t 




radius of hole in insulator layer 


L 




length of the propellant grain 


L 




typical length equal to lm 


s b 




total area of the burning surface 


F N 




normal thrust 


F h 




additional thrust produced by hole gas flow 


r b 




burning rate of solid propellant 


n 




exponent for burning rate of the propellant 


Pp 




density of the solid propellant 


H P 




combustion heat of the solid propellant 


Q 




heat flow from the gas to the walls of the hole 


s 




cross-section of the combustion chamber 


ftr 




surface friction force 
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T m = melting temperature point 

T„ = critical temperature of the nozzle ablation 

T c = temperature of metal case far from hole 

c m = specific heat of case metal 

q ins = latent heat of insulator ablation 

q m = latent heat of metal melting 

c„ = specific heat of nozzle material 

q„ = specific heat of ablation of insulator layer 

p m = density of case metal 

p„ = density of nozzle material 

k = the thermal conductivity 

/j = dynamical viscosity of hot gas 

Pr = the Prandtl number, PT=^C p /k 

h = subscript for gas parameters in the hole 

N = subscript for parameters in normal regime 

= subscript for stagnation values of gas parameters 

1. Introduction 

New generation of the heavy lift vehicles for space 
exploration, including ORION mission, requires 
development of the new IVHM systems with the 
overarching goal of safer and more reliable flights. The 
difficulties in the development of an on-board IVHM for 
SRBs are the facts that: (i) internal hydrodynamics of SRBs 
is highly nonlinear, (ii) there is a number of failure modes 
that may lead to abrupt changes of SRBs parameters, (iii) 
the number and types of available sensors are severely 
limited, and (iv) the recovery time is typically a few 
seconds. These difficulties dictate the model based approach 
to the IVHM of SRBs that minimizes the number of 
"misses" and "false alarms" [4] by utilizing deep 
understanding of the physical processes underlying the 
nominal and off-nominal regimes of SRBs. 

Indeed, theoretical analysis shows that many fault modes 
leading to the SRBs failures [5]-[7], including combustion 
instabilities [8]-[10], bore choking [11]-[13], and case 
breach [5], [14], have unique dynamical features in the time- 



traces of pressure and thrust. The corresponding expert 
knowledge can be incorporated into on-board FD&P within 
the novel Bayesian inferential framework [l]-[4] allowing 
for faster and more reliable identification of the nominal and 
off-nominal regimes of SRB operation in real time. 

However, the development of such an inferential framework 
in practice is a highly nontrivial task since the internal 
dynamics of gas flow in the SRBs is essentially nonlinear 
and results from the interplay of a large number of complex 
nonlinear dynamical phenomena on the propellant and 
insulator and metal surfaces, in the gas flow, and in the 
nozzle. On-board FD&P, on the other hand, can only 
incorporate low-dimensional performance models (LDPM) 
of the pressure and thrust dynamics of SRBs. The derivation 
of the corresponding LDPMs and its validation in high- 
fidelity simulations and experimental rocket firing tests are 
important ingredients of the development of the FD&P. 

A very important issue of the development of on-board 
IVHM system for SRBs is its ability to accommodate the 
abrupt changes of the model parameters in various nonlinear 
dynamical off-nominal regimes. The difficulties in building 
model-based hybrid probabilistic algorithms for on-board 
IVHM stem from the fact that at present there are no general 
inferential frameworks for nonlinear dynamical systems 
with mixed binary and continuous parameters. 

In the present paper we address both problems in the context 
of derivation and testing of a LDPM for on-board FD&P. 
Specifically, we introduce a high-fidelity model of the case 
breach fault, derive a corresponding LDPM, and compare 
the results of the simulations of both models. The LDPM is 
then applied to analyze two situations motivated by recent 
ground firing tests. In the first example, we assume that the 
hole in the metal case appears suddenly (perhaps due to an 
interaction with an external object), the hole in the metal is 
larger than the hole in the insulator and the dynamics of the 
fault is determined by the dynamics of the ablation of the 
insulator. We show that, in this case, FD&P can be 
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Figure 1 (left) Geometry of the model surfaces in FLUENT. Note that the hole wall, propellant surface wall, and the 
nozzle wall are deforming. The dashed blue lines on the left and right hand side show boundaries of the ambient 
regions, (right) Velocity distribution generated by the FLUENT model for t = 0.34 sec. 
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Figure 2 Axial velocity (left) and pressure (right) profiles generated by the FLUENT model for t = 0.05 sec (black 
solid line) as compared to the analytical solutions (red dashed lines) given by the Eqs (7)-(8). 



performed using a self-consistent iterative algorithm that 
avoids numerical integration of the LDPM. It is assumed in 
this case that the geometry of the fault is simple, so that the 
dynamics of the fault is governed by equations for the 
ablation of the nozzle and insulator, and melting and 
burning of the hole in the metal case based on Bartz' 
approximations [20] -[23] . 

In the second example we consider the situation when the 
hole in the metal case develops more slowly perhaps due to 
the crack in the insulator. In this case the dynamics is 
governed by the melting of the metal. The geometry of the 
fault in this example is more complicated and is motivated 
by the recent ground tests. In this example the fault 
dynamics deviates substantially from the one given by the 
Bartz' approximation. Accordingly a more general 
algorithm for FD&P is introduced. 

We show that the introduced algorithms can accommodate 
abrupt changes in the model parameters of essentially 
nonlinear dynamics of the gas flow in SRBs. The 
applicability of the results to different SRB fault modes is 
discussed. 

2. Case breach model 

In this section we develop both axisymmetric FLUENT [25] 
HMF of the case breach fault and a set of PDEs that 
describe this fault in ID approximation. 

To model a case breach fault (see the sketch in Figure 1) we 
assume that the gas leak through the hole in the case is small 
compared to the gas flow through the rocket nozzle. An 
analysis of the fault caused by a small gas leak is the most 
important one from the point of view of the development of 
prognostics and diagnostic for a SRB. An axisymmetric 
HFM of this fault is built in FLUENT. The geometry of the 
model surfaces is shown in Figure 1 . This model takes into 
account the nozzle ablation and the melting of the metal in 
the hole through the rocket case. Even more importantly, it 
models the dynamics of the burning area as a given function 



of the burn web distance. These features allow us to 
reproduce accurately experimentally observed time -traces of 
the rocket and gas flow parameters. An example of the 
velocity distribution generated by this model is shown in 
Figure 1. 

The dynamics of the distributions of the gas flow parameters 
obtained with the FLUENT model can now be applied to 
validate low-dimensional performance model (LDPM) of 
SRBs. To do so we notice that in axisymmetric case one can 
use one-dimensional approximation for both gas flow in the 
combustion chamber and gas flow in the hole 1 . To describe 
gas flow in the combustion chamber in ID approximation 
we can use a set of equations for the dynamics of mass, 
momentum, and the energy of the gas, introduced in our 
earlier work [3], [4]. 



8 ( (Spu) = -8 X {spu 2 ^ - Sd x p + A 2 ^ 2 (t), 
d t ( S P e t ) = ~ d x ( s P h t u ) + H Pp r b l + (0- 



(1) 



Here d^dldt, d x =dldx, e,=c v T+u /2, h,=c P T+u /2 are the total 
energy and total enthalpy of the gas flow in the combustion 
chamber, £(t) are white Gaussian noises with intensities A t . 
The boundary conditions for Eqs (1) assume the continuity 
conditions at the nozzle inlet and the stagnation values for 
the gas parameters at the point where the Mach number is 
zero in the combustion chamber. Eqs. (1) represent a 
modification of the well-known model of Salita [15]. Note 
also that in general case the gas dynamics is stochastic and 
the noise terms A^t) should be included into the equations. 

We now extend this mode by coupling the gas dynamics in 
the combustion chamber to the gas flow in the hole. The 
corresponding set of PDEs 



1 Note that the same approximation is valid for arbitrary geometry of 
the fault in the limit small gas leak when the perturbation of the gas flow 
in the combustion chamber is small. 
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8 t{ s hPh) = - 8 x{ s hPh u h\ 

d t [ s hPh u h ) = ~ d x [ s hPh u h ) - s h d xPh - f fr l h ' ( 2 ) 



t\ s h p h e t 



s h.Ph u h? l t 



h l h> 



resembles Eqs. (1). The important difference, however, is 
that we neglect mass flow from the walls of the hole. 
Instead Eqs. (2) include the term that describes the heat 
flow from the gas to the hole walls. The boundary 
conditions for this set of equations assume ambient 
conditions at the hole outlet and the continuity equation for 
the gas flow in the hole coupled to the sonic condition at the 
hole throat. The value of Q h is presented in the Appendix. 
The dynamics of the gas flow in the nozzle is described by a 
set of equations similar to (2) and can be obtained from this 
set by substituting subscript "n" for subscript "h". 

The sets of PDEs for the gas dynamics in the combustion 
chamber (1), hole (2), and nozzle area are coupled with the 
equation for the steady burning of the propellant 



dR 
dt 



f \ 
P_ 

\Pc 

and equations for the nozzle ablation 
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and melting of the metal walls of the hole in the rocket case 
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Velocities of the melting v m and ablation v„ fronts are given 
in the Appendix. 

It is further assumed that throughout the combustion 
chamber and in the hole the following equation of a perfect 
gas holds 



P=( Cp -C r )T=P± 
P Po 
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(6) 



In a quasi-steady regime Eqs (l)-(6) can be integrated 
analytically [4], [14] (cf. with numerical solution in [16]) 
leading to the following set of equations in the combustion 
chamber 
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In the nozzle and the hole areas the analytical solution can 

be written in a well-known (see e.g. [17]) form 

P = Po M fact> P = Po M fact> 



T = T n M f t , 
fact' 
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where M fact = (l-(f-l)M 2 /2). The results of the 

analytical solution are in good agreement with the results of 
the CFD simulations in FLUENT as shown in Figure 2 for 
time t = 0.05 sec. Similar results are obtained for later times, 
demonstrating a good agreement between theory and 
simulations of time variation of quasi-steady axial 
distributions of the flow parameters. The good agreement 
between the analytical results and the time variation of the 
axial distributions of the gas flow parameters in FLUENT is 
the first step in the validation of the low-dimensional 
performance model of the fault derived in the next Section. 

3. Low-dimensional performance model 

It follows from the results of analytical calculations and 

simulations of the FLUENT HFM that M 2 =v 2 I c\ = 1 

is small everywhere in the combustion chamber. 
Furthermore, the equilibration of the gas flow parameters in 
the combustion chamber occurs on the time scale of the 
order of t = Lie « 6 msec. As a result, the distribution of the 
flow parameters follows adiabatically the changes in the 
rocket geometry induced by the burning of the propellant 
surface, nozzle ablation and metal melting in the hole 
through the case. It is, therefore, possible to derive the 
LDPM of the case breach fault by integrating Eqs (1), (2) 
(see [3], [4], [14] for further details) and combining the 
result of integration with Eqs (3)-(6) to obtain the following 
set of six ODEs 
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Here the following dimensionless variables are used 
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where subscript m refers to maximum reference values of 
the pressure and density. Important novel features of the 
case breach LDPM derived above are the following: 

(1) The burning area of the propellant is calculated 

using a given design curve S^= f (R) that relates 

to the burn web distance R. 

(2) The dynamics of the case breach fault is taken into 
account by including the last equation for melting the 

metal walls of the hole through the rocket metal case. 

2 2 

(3) The effective nozzle area s et = Tir h&w t is 
introduced in the first two equations that take into 
account the effect of the fault dynamics on dynamics 
of the flow parameters at the stagnation point. 

(4) The dynamics of the volume of the combustion 
chamber is taken into account in the equation 4 of the 
set of Eqs. (9). 

(5) The dynamics of the nozzle ablation is taken into 
account in the last equation of the set. 

In the general case the set of Eqs (9) may also include an 
equation for the dynamics of the hole in the insulator layer 

f \~Pi 



1-A 



(11) 



These features allow us to model realistic changes in the 
rocket geometry with arbitrary thrust curves. 

In the following sections we discuss probabilistic algorithms 
for an on-board SRB FD&P based on the model introduced 
above. These algorithms should satisfy a few very strict 
requirements mentioned in the introduction. In particular, 
these algorithms should accommodate very short time- 



windows for the diagnostics and prognostics (of the order of 
a few seconds) and abrupt changes of the model parameters 
reflecting sudden transitions from nominal to off-nominal 
regimes. As the main input, they use time-traces of 
stagnation pressure and thrust. These algorithms should 
allow for discrimination between multiple fault modes of 
SRBs. 

In particular, we consider two algorithms. The first 
algorithm is based on the self-consistent algorithm that 
solves the LDPM in a quasi-adiabatic approximation, when 
the pressure and density follow adiabatically dynamics of 
the propellant burning, melting of the metal case, and nozzle 
and insulator ablation. The second algorithm is based on the 
integration of the system of stochastic differential equations 
of the LDPM. 

4. Inference of the model parameters 

Note that effect of the case breach fault on the dynamics of 
the internal gas flow in SRBs is reduced to the effective 
modification of the nozzle throat area as explained above. It 
is, therefore, possible to infer SRB parameters using a 
Bayesian framework introduced in our earlier work [l]-[4], 
[14], [30] for an analysis of the overpressure faults due to 
the changes of the nozzle throat area. In particular, it was 
shown [4], [14], [30] that this algorithm can accommodate 
sudden changes of the model parameters and, therefore, is 
suitable for developing of the hybrid probabilistic IVHM of 
SRBs. 

Here we briefly reproduce earlier the results related to the 
analysis of the abrupt changes of the model parameters. The 
dynamics of the LDPM (9) can be in general presented in 
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Figure 3 Estimation of the value of the parameter - 
CQrSf/faL) before (left curve) and after (right curve) 

the step-wise change of the nozzle throat area. The 
dashed line shows the actual value of the parameter. 
The solid lines show the PDF of the parameter 
estimation with T=0.1 sec, Af=h=0.001 sec, ^=500. The 
inset shows the sudden change in the nozzle throat 
geometry. 
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the form of Euler approximation of the set of ODEs on a 
discrete time lattice {t k = hk; k = 0, 1, . . ., K} with time 
constant h 

Jhz k , (12) 



X M= X k+¥(x k \c) + 

1 r'i, +h * 
where z k =—=^ E,(t)dt , x k 



X k + X k+1 



{p,p,R,V,r h ,r t ,r^ is Z, -dimensional state of the system (9), 
(11), <r is in our case a diagonal noise matrix with two first 
non-zero elements aj and a 2 , f is a vector field representing 
the rhs of this system, and c are parameters of the model. 
Given prior distribution for the unknown model parameters 
in a form of Gaussian distribution we can apply our theory 
of Bayesian inference of dynamical systems [l]-[4] to 
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Here the vector field is parameterized in special form [1] 
f(x\ c) = U (x)c , where U(x) is a block-matrix with N 

blocks of the form I(/) n {x(t k y) , / is LxL unit matrix, and 

v ™( x ) = X^ — • 

„=i ox„ 
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Table 1 Results of the parameter estimation of the model (9) 
in the absence of the nozzle ablation. The total time of the 
measurements in this test was 7M).5 sec, the sampling time 
was A?=0.001 sec, and the number of measured points was 
N=500. 

The application of the algorithm (13)-(16) to the problem of 
inference of the abrupt changes of the gas flow parameters 
in the SRB due to the changes in the rocket geometry are 
presented in Figure 3 and Table 1. It can be seen from the 
figure and table that the algorithm can (at least in principle) 
infer sudden changes in the gas flow with time resolution of 
the order of 0.1 sec with relative error less than 1%. 

Below it is assumed that the algorithm was applied to learn 
the SRB parameters in the nominal regime and multiple 
ground and flying tests. Therefore, nominal SRB parameters 
are considered to be known and we are mainly concerned 
with diagnostics and prognostics of the case breach fault 
parameters. 

5. Example I: Self-consistent iterative 

ALGORITHM 

In this example we analyze the following problem. A hole 
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Figure 4 (left) Results of the self-consistent calculations of the rocket parameters in the off-nominal regime using 
iterative algorithm A. Absolute values of pressure for four different initial values of the hole in the insulator: 0.5, 1.0, 
1.5 and 2.0 mm are shown by the black, blue, red, and cyan solid lines respectively. The nominal pressure is shown by 
the dashed black line, (right) Iterations of the effective hole radius in the metal case. th approximation is shown by 
the red solid line. Five first approximations shown by red dashed lines are indicated by arrows. Final radius of the 
hole in the metal case is shown by black line. th approximation for the hole in the insulator is shown by dashed blue 
line. Final radius of the hole in the insulator is shown by the black dashed line. 
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t, sec t, sec 

Figure 5 (left) Fault-induced thrust (black solid line) is shown in comparison with nominal SRB thrust (blue solid line) 
and off-nominal SRB thrust (dashed blue line). Initial radius of the hole in the insulator is 0.75mm. (right) Pressure 
(black line) and temperature (blue line) in the metal hole through the case determined by the iteration algorithm B. 



through the metal case and insulator occurs suddenly at the 
initial time of the fault t . We have to infer and predict the 
dynamics of the growth of the holes in the insulator layer 
and in the metal case, as well as the fault-induced side 
thrust, and changes of the SRB thrust in the off-nominal 
regime. As an input, we use time-traces of the stagnation 
pressure in the nominal regime and nominal values of the 
SRB parameters. In particular, it is assumed that the ablation 
parameters for the nozzle and insulator materials and the 
melting parameters for the metal case are known. In this 
case we assume that the hole radius in the metal case is 
larger than the hole radius in the insulator (i.e. the velocity 
of the ablation of the insulator material is smaller than the 
velocity of the melting front), accordingly the fault 
dynamics is determined by the insulator ablation of the 
insulator (see details in [3], [14]). 

To solve this problem we introduce a prognostics algorithm 
of the fault dynamics based on a self-consistent iterative 
algorithm that avoids numerical solution of the LDPM. We 
notice that with the limit of steady burning, the equations in 
(9), (11) can be integrated analytically. Because the hole 
throat is determined by the radius of the hole in the insulator 
layer, at this stage we can omit the equation for the hole 
radius in metal case. The resulting set of equations has the 
form 



P(t) = P N (t) 



S bef^ S tN^ 
S bN s et^ 
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l-n 



R-a\ p (z)dr, 



(17) 
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Here S et (t)=S tN (t) + As t (t)+S h (t)is an effective 

nozzle throat area where the 1 st term corresponds to the 
nominal regime, the 2 nd corresponds to the deviation of the 
nozzle throat area from the nominal regime due to the fault, 
and the 3 rd term corresponds to the area of the hole throat in 
the rocket case. Similarly, we define the effective burning 

area S bet (f) = S bN (f) + AS b (t) as a sum of the burning 

area in the nominal regime and a term that describes the 
deviation of the burning area from the nominal regime due 
to the fault. Using Eqs (17) the following iterative algorithm 
Al can be introduced: 

(1) Set initial values of the corrections to the nozzle and 
burning area to zero As,(t) = and ASb(t)=0. Set values 
of areas of the holes in the insulator and in the metal to 
constant initial values si,(t) = n-r u 2 and s met (t) = n-r mel0 2 . 

(2) Update time -trace of the pressure using 1 st Eq. in (17) 

(3) Update burn web distance R, radius of the hole in the 
insulator r h , and nozzle throat radius r t , 

(4) Repeat from the step (2) until convergence is reached. 

The results of the application of this self-consistent 
algorithm to the prognostics of the case breach fault 
parameters are shown in Figure 4 (left). Once quasi-steady 
pressure and the dynamics of growth of the hole in the 
insulator are predicted in the off-nominal regime one can 
determine the dynamics of the hole growth in the metal case 
and the dynamics of the fault-induced side thrust. To do so, 
we use the following self-consistent iterative algorithm B 
for t>t that takes into account the assumption that the 
velocity of the melting front is larger than the velocity of 
ablation in the insulator 

(1) Set th approximation r (0> met(l (t) for the hole radius in the 
metal to r met0 . 

(2) Construct 1 st approximation 
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The deep physical meaning of the iterative procedure 
introduced above is that the ablation of the nozzle and the 
case breach fault develop in a self-consistent manner. 
Indeed, the increase of the leak cross section, due to 
insulator ablation under the action of the hot gas flow, leads 
to decreased pressure and hence a decreased burning rate. 
This, in turn, decreases the hot gas flow through the hole 
and the ablation rate. In this way a quasi-stationary regime 
of burning and ablation is developed. The parameters of 
burning in this regime can be found in a self-consistent way 
using an iterative algorithm, without integration of the full 
system of differential equations of motion. 



dr 



6. Example II: Case Breach Fault with 
Nontrivial Geometry 



(5) 



(6) 
(7) 



Use Eqs (8) to find velocity u h (t)=M oh (t)-c , mass flow 
(pu) h (t), temperature T h (t), and pressure p h {t) in the 
metal hole; 

Repeat from the step (2) until convergence is reached. 
Calculate fault-induced side thrust 



(pu) 



h u h,ex + (Pfi,ex 



Pa) 



met 



A similar algorithm is used to find the ablation of the nozzle 
and SRB thrust in the off-nominal regime in the case breach 
fault. The results of these calculations are shown in Figure 4 
(right) and Figure 5. We note that the inference of the fault 
parameters is trivially achieved using the same iterative 
algorithm with the only exception being that the time-trace 
for pressure in the off-nominal regime is given by the 
measurements. Accordingly the first equation in the set (17) 
and the 2 nd step in the iterative algorithm A are not needed. 
Note also that the important part of the introduced above 
algorithms is an assumption that the design curve Sy= f[R) 
representing the relation between the burning area S A and 
burn web distance R is know and remains invariant 
characteristics of the SRB in the off-nominal regime of the 
case breach fault. 



In this Section we consider a situation, motivated by recent 
ground tests, in which fault occurs in a closure holding a 
pressure sensor. The disk has a form of a metal disk of 
radius Rj and is covered by insulator material. We assume 
that a small leak arises on a boundary of the disk at moment 
t=t (see Figure 6). Subsequently the small hole will grow 
due to surface melting and flame erosion (burning) of the 
metal under the action of the hot gas flow effluent from the 
combustion chamber. The problem of the FD&P is to detect 
fault, infer a model of the fault dynamics, and to predict 
pressure and thrust dynamics ahead of time in the presence 
of the fault. At the input of the algorithm we have time- 
traces of pressure and thrust, and known nominal parameters 
of the SRB and ablation of the nozzle. The difficulty of this 
problem is related to the fact that in the presence of the 
complicated geometry of the fault it is not sufficient to 
know the melting and burning parameters of the metal case. 
Indeed, in this case the last equation in (9) has to be 
substituted with the following set of equations 
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1200 





Figure 7 (left) Pressure time-trace p(f) in case breach fault-induced off-nominal regime (blue line) as compared to p(t) 
in the nominal regime (black line), (right) The result of the reconstruction of the time-trace of the hole area s h (t) (blue 
line) using algorithm C as compared to the actual dynamics of the fault (black line). 



Qc+Qm+O, 



\jm+ C m { T mel- T m o)\pi 

Qc + Qr 



+ V ft" 



(18) 



n - {l- x^jarccos 




Here we take into account the geometry of the fault and the 
heating of the disk due to the convective and radiation flows 

with the values of Q R and Q c given by equations (29) and 

(32) and V/j is a constant rate of erosive burning (see 
Appendix). Note that to find ablation of the nozzle exit we 
can use time-traces of p(f) and r,(t) obtained by integration 
of 3 equation in set of Eqs. (9) and the equation for the 
ablation of the nozzle exit 



order of 0.1 in/sec, which is frequently observed in 
experiments (see e.g. [29]). The reason is that there exists a 
delay time (related to the heating of the metal) before a 
steady regime of melting can be established. As a 
consequence the typical time-trace of the pressure in the off- 
nominal regime may have the form shown in the Figure 7 
(left) by the blue line. In this case the hole in the case 
appears at to = 1 sec and grows with a velocity of 
approximately 0. 16 in/sec. 

To infer parameters of the fault we introduce the following 
algorithm C 



&= v 




-fit 



T-T 
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(19) 



(1) 

(2) 
(3) 
(4) 



a J 



where T„ is given by the solution of the last two equations 
in (8) for a given S,=nr„ 2 . 

Note that the relation between hole radius r h and hole area s h 
was also changed to take into account the complex geometry 
of the fault. Accordingly, the fault dynamics becomes 
highly nonlinear. This is true even in the case of metal 
melting with constant so-called erosion velocity v er of the 



Given the pressure time-traces p(i), integrate 
equations for the nozzle ablation (equation 5 in (9) and 
equation (19)) to obtain s,(t) and s ex {f) in the off- 
nominal regime. 

Use equation 3 in (9) and p(t) to obtain time-trace of 
the web distance R(t) 

Use R(t) and the design curve Sb= f(R) to find time- 
trace of the burning area Sb(t) 

Use Sb(t) and quasi-stationary solution of the equation 
for the pressure to find a time trace of the effective 
area of the nozzle throat 



"et 



(t) 



Pit) 



Poo ' Poo 



J 



aCoPpPpT 

r 



(5) Find time -trace of the hole throat area Sh{t)=s et {t)-s t {f) 

(6) Find nozzle thrust and fault-induced thrust. 
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Figure 9 (left) Actual fault dynamics is shown by the blue solid line. The time interval elapsed from the case 
breach fault is shown by the black circles. The predicted trajectories of the fault are shown by thing green solid 
lines. The time moment of the prediction is indicated by a vertical red line, (right) The distribution of the 
predicted value of the fault at the future time t = 8 sec is shown by solid blue, black, and green lines in 
comparison with the true value of the fault indicated by the vertical red dashed line. The time elapsed from the 
case breach fault used for predictions is shown in the figure. 



The result of the application of this algorithm to the 
analysis of the time-trace of the pressure is shown in Figure 
7 (right) and Figure 8. It can be seen from the figures that 
the algorithm C can infer accurately the fault dynamics. We 
now apply these results to develop an algorithm of 
prognostics. 

7. Prognostics algorithm for complex 
geometry 

First we note that prognostics algorithms are much more 
sensitive towards the errors of parameter estimations than 
the diagnostics algorithms. For example a simple 
polynomial fit of the first 2.5 sec of the inferred time -trace 
of si,{t) shown in Figure 7 (left) is not sufficient to predict 
fault dynamics. The corresponding predictions will quickly 
diverge from the actual value of Sh(t) » 2 sec later. 
Therefore, more sophisticated algorithms are required to 
predict nonlinear fault dynamics. Such algorithms can be 
developed using the LDPM (9) and Bayesian dynamical 
inference discussed briefly in Sec. 4. 

To apply the dynamical inference algorithm to the 
prognostics of the case breach fault we notice that metal 
erosion of the hole walls during the case breach fault is 
mainly determined by the constant rate Vj (see e.g. [29], 
model (18), and the results of simulations shown in Figure 7 
(right)). And if the fault geometry is simple, which 
corresponds to the hole in the side wall, the si,{t) oc t, and the 
predictions will be considerably simplified. If, however, the 
hole occurs at one of the joins, which is very likely the case 
in practice, the time evolution of the S/,(0 will be mostly 
determined by the geometry of the joins and the predictions 



will be mainly based on the earlier detection of the 
catastrophic trends in the Sf,(t). 

The catastrophic trends in the s/,(t) are usually related to the 
high powers f of t, because these terms determine fast 
deviations of the dynamics from the nominal regime. In this 
context it is important to verify convergence of the method 
that includes higher order terms into the set of the base 
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Figure 8 Nominal pressure (black line) as compared 
with the off-nominal pressure of the nozzle (blue 
dashed line) and of the hole (blue dashed lines at the 
bottom). Inferred values of the nozzle and hole thrusts 
are shown by red dashed lines. 
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functions of s h (i). 

To model such situation we will consider practically 
important case when only pressure signals p(t) (and perhaps 
temperature T(t)) are given by the measurements. Then the 
model (9) can be reduced to the following form 

^- C -^pjT+^{yp p -p)p n +jDtM (20) 

Vr b V v ' 

Here T ' = pip » 1 for dimensionless pressure and density 
introduced in (10), and the time traces R{t), F(t), V(t), and 
r t (t) are completely determined by the pressure and inferred 
law of the nozzle ablation 



R(t) = \[p n (V)dt\ F(t) = f(R(t)), 
V(t) = V +\l (,) f(R)dR, s t (t) = nrfitX 



(21) 



r, (t) = 



'to ^ 



up 



and, therefore, are known. 

We choose s h (t) = s et (t) - s t (t) in the form 



(22) 



m—\ 



and apply dynamical inference algorithm introduce in Sec. 4 
to infer parameters a, and D in (20). All other parameters in 
(20) are know from the analysis of the nominal regime. 
The ID problem (20), (22) can be solved by applying Eqs 
(13)-(16), where U(x„) has the form 

U = {l,t,t 2 ,t\t 4 ,t 5 }. (23) 

The results of the predictions are shown in the Figure 9. In 
the left figure we show the time elapsed from the fault (t = 2 
sec) by black circles. The pressure time -traces measured 
during the first 2 sec are then used to infer parameters a, and 
D. Next, parameters a, are used to predict the future time 
evolution of the fault using (22). The distribution of the 
predicted trajectories obtained in this way is shown by the 
blue line in the Figure 9 (right). When new data are 
measured on the interval 2.5 sec and 3 sec the procedure is 
repeated to obtain the updated distribution of the predicted 
values of the fault as shown in the Figure 9 (right) by black 
and green lines. It can be seen from the figure that the 
method converges during first 3 sec of the fault. 



number of algorithms are derived to infer fault parameters 
and to predict fault dynamics. The algorithms accommodate 
abrupt changes in the model parameters and can be used to 
develop hybrid probabilistic on-board SRB IVHM. The 
performance of the algorithms was tested using analysis of 
the experimental time-series data. It is shown that the 
algorithm can be successfully applied for the prognostics of 
the case breach fault. 

The developed methods and algorithms can be used to 
analyze other SBR' faults, including overpressure and 
breakage of the case induced by nozzle blocking, bore 
choking and grain deformation. The bore choking 
phenomenon is an almost radial deformation (bulge) of the 
propellant near booster joint segments. This phenomenon 
can cause choking of the exhaust gas flow and increase the 
burning surface which can lead to critical overpressure in 
the combustion chamber. Development testing has shown 
that this fault was observed, for example, in primary 
construction of the Titan IV (see the report [6], [7], [31]). 
The bore choking, and also cracks and voids in the solid 
propellant, can result in local burning of the booster case 
and also in abrupt breaking off of large pieces of the 
propellant. These pieces can stick in the nozzle throat and 
block the exhaust gas flow (nozzle blocking fault). In all 
these cases the fault dynamics is governed by the changes in 
the burning area and/or effective nozzle area. Therefore the 
LDPM introduced in this work can be efficiently applied to 
analyze of these faults and develop on-board SRB IVHM. 

Appendix : 

The goal of the appendix is to study ablation of the nozzle 
surface and also melting and burning of the metal surface of 
a leak hole. These surfaces heat under the action of a hot gas 
stream Q flowing through the hole or nozzle. The surfaces 

start to corrode when the surface temperature T s reaches the 
ablation temperature of the insulator material T a u or the 
melting point of the metal T md , respectively. Let us analyze 
the growth of T s on the surface of a hole of a radius R h in a 
metal disk or nozzle throat of a radius R N . Temperature 
distribution inside the material is described by 



dT 
dt 



( a2 



K 



V 



d z T 



1 dT_ 
r dr 



(24) 



8. Summary and discussion 

In this paper we report the work in progress on the 
development of hybrid probabilistic algorithms for SRBs 
on-board IVHM system. Two experimentally motivated 
situations are considered. In the first situation case breach 
fault occurs suddenly in the metal case with simple 
geometry. In the second situation the fault develops slowly 
in the forward closure with complex geometry. A high- 
fidelity model of the fault dynamics is developed in 
FLUENT. The LDPM of the off-nominal regime induced by 
the case breach is derived using the high-fidelity model. A 



where r is a coordinate normal to the surface at r>i? hN . It is 
well known that a homogeneous solution of this equation is 

2 

where y 1 = — (25) 

Cp 



T(x, t) = A{t) exp 



4y 2 t 



i. e. thermo-diffusion length is equal to: 

l D = l^Ktl cp (26) 
Here the parameters C, p, k correspond to a metal or an 
insulator material of a nozzle. The flow of the hot gas Q for 
the time interval dt will heat the surface layer of an area dS 
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and a thickness of l D . For example, the heat balance for the 
metal surface can be written as: 



Q c = 0.023? 7r C p 



yp 



^ •^ 2i^ ^-o. 2 



v rc oy 



(T -T s ) 



(32) 



QdSdt-cpdSd[(T s -T m0 )l D j (27) where 7* = 27 / (l + 7) . The increase of the surface 



V M J 
The 



where T m0 is an initial temperature of the metal disk. It 
follows from Eqs. (26) and (27) 



dT c 



( T s- T mo) 



(28) 



Q 

dt It 2 N C'/)A7 

Here we assume Q = Q h =Q C +Q R where Q R and Q c are 

radiation and convective heat flows from the hot gas to the 
metal surface. Radiation heat flow is given by 



where A 



Q R =a[l- SW (-AP)](T 4 -T* ) (29) 

0.001 + 4 x 10" 4 x%^l!/ 14.69 /psi is the 

emissivity of the hot gas, %AL is the percentage of 
aluminum in a solid propellant, and 
-8 2 4 

a = 5.67x10 Wlm K is the Stefan-Boltzmann 
constant. For realistic parameters XP « 1 , therefore 

Q R ; gXp[t 4 -T^q) . The convective heat flow 
Q c =hg{T-T s ). When the surface temperature T s 

reaches to T mel , the value of T s is equal to T } inasmuch 

as the melted metal layer rushes out from the hole almost 
immediately for a time Klmsec [14]. The heat transfer 
coefficient h g can be written as [20] -[23]. 



h 



(30) 



: 0.5JC p r] r C f 

Where the coefficient rj r describes an effect of the 

roughness. The value of the coefficient Q is determined by 
Reynolds number: Re = pux I /u . Re » 10 6 for our 

parameters and p=1000psi. The critical Reynolds number 
which determines a transition between laminar and turbulent 
regimes is known to vary from approximately 10 5 to 3xl0 6 , 
depending on the surface roughness and the turbulence level 
of the free stream [20]-[21]. The turbulence regime probably 
occurs even for p ; 100 psi due to growth of the hole radius 
and the roughness of the metal surface. Q can be 
approximated over a fairly wide range of Reynolds numbers 
by[20]-[23] 



Cj = 0.046 



JD 



aa^d -°'2 n "0.67 

0.046 Re Pr 



(31) 



Here the Prandtl number Pr = juC p I k ; 1 [19]-[21]. It 
follows from Eqs. (30) and (31) 



temperature is given by Eq. (28) where Q c and Q R are 

given by Eqs. (32) and (29). Eq. (28) determines a time t=to 
when the surface temperature T s reaches T m . Our 
calculations show that to<0.2sec for a steel and t «0.1sec 
for a typical nozzle material. 

The melting starts at time t>t when T s =T md . The velocity of 
melting can be found from the heat balance equation: 

QdSdt = [q + c(T-T )]pdSv f dl D (33) 

where T is an initial temperature. It follows that the 
velocity of melting front for time t>t is given by 



*=v, = r a+Q « n 

me ' f [q m + C m (T mel -T )]p m 



(34) 



The metal surface can react with oxidizing agents of the 
combustion products (hot gas). Obviously, velocity of 
propagation of metal burning front vj^ has to become 

saturated at high pressures and/or intense surface flows of 
oxidizing agents. Experiments show that 
vjfo ; 0.76cm / sec for iron and carbon steels [21], [22]. 

The same value of vr^ was reported in the study of ATK 

Thiokon Inc. team [29]. We notice that the gas temperature 
T=T* inside the hole is greater than temperature of the 

metal surface, therefore a heat flow Q h induced by the 

burning spends for the metal melting. Taking into account 

the burning of the metal surface we can written v ^ in the 

following form 

& = v Qc+Qr+Q, 

f ~ll m +C m (T mel -T m0 )]p m (35) 

= Q^Qm +v 

[l m +C m (T mel -T m0 )]p m fi 

Analogously Eq. (34), the velocity of insulator material 
ablation of nozzle t>t is given by 

Qc+Qr 



(36) 



[q in + c in (T ahl -T )] Pin 

Estimations show that Q c is much greater than bothg^ 

Thus, we can use Eqs. (4) and (5) where v n and v m are 
approximately equal to 
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